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Abstract 

Monte Carlo simulations and an analytical approach within the framework 
of a semiclassical model are presented which permit the determination of 
Coulomb blockade and single electron charging effects for multiple tunnel 
junctions coupled in series. The Coulomb gap in the I(V) curves can be 
expressed as a simple function of the capacitances in the series. Furthermore, 
the magnitude of the differential conductivity at current onset is calculated 
in terms of the model. The results are discussed with respect to the number 
of junctions. 
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Since the development of quantum mechanics, electron tunneling has been widely inves- 
tigated experimentally and extensively discussed in theory. |||4[] The first experiments 
were performed on metal-insulator-metal (MIM) sandwiches, [[I] which permitted the study 
of tunneling phenomena and the direct correlation of such effects with the thickness of the 
potential barrier. There are two main theoretical approaches commonly used to describe 
tunneling phenomena. The first entails solving the time-independent Schrodinger equation 
and matching the wave functions. [|J However, an exact matching of the wave functions can 
only be achieved in a one- dimensional approximation. The second method is known as the 
transfer-Hamiltonian approach, Q which generally involves a time-dependent Schrodinger 
equation in the second quantization formalism. In this approach, tunneling of the electrons 
between the two electrodes is described using a transfer Hamiltonian in addition to the 
Hamiltonians of the unperturbed electrodes. 

For small-scaled systems a particular tunneling phenomenon referred to as single electron 
tunneling (SET) || can occur. SET may be observed, for example, in a system of two 
electrodes separated by a thin insulating layer, in which a third electrode (e.g., a small 
metal particle) is embedded. The current through the system is controlled by a single 
electron if k B T <C e 2 / [2(Ci + C 2 )], with C\ and C 2 being the capacitances of the two 
junctions. A tunneling event therefore changes the total charge of the center electrode by 
±e, depending on the direction of tunneling, and thus the electrostatic energy of the system. 
At zero temperature, tunneling is completely prohibited for \V\ < e/ [2(C\ + C2)} = AV^ ap /2 
(Coulomb blockade), and the transport process is dominated by charging effects. This 
means, in particular, that an increase in voltage applied to the double junction system 
leads to incremental charging, which might manifest itself as steps in the I{V) characteristic 
(Coulomb staircase). || Such I(V) characteristics are usually measured at temperatures 
below 4 K [0 for certain relations of RiCi, such as -R1C1 <C R 2 C 2 (Ri and Cj are the tunnel 
resistance and capacitance of the i-th junction). The theory developed so far only considers 
this low-temperature limit. 

One important point in understanding SET is the extension of the theory to multiple 
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tunnel junctions as well as to room temperature. Increasing the temperature to room tem- 
perature demands capacitances as small as C — 0.1 aF, since only then does the electrostatic 
energy of charging by a single electron, e 2 /2C, exceed the thermal energy /c#T \ Nowadays, 
such capacitances can be realized using small particles |7| or clusters || with dimensions in 
the 1 to 10 nm range. 

Here, a theoretical study of I(V) characteristics of one- dimensional chains of small ca- 
pacitors coupled in series is presented. Such arrangements can be realized, for example, 
in granular films of small metal islands or with clusters. An analytical approach for zero 
temperature in the framework of a semiclassical model is used to explain several features 
arising in Monte Carlo simulations at finite temperatures. It is shown that the width of the 
Coulomb gap is directly related to the number of capacitors in the series. An expression 
is given for calculating the voltages at which steps in the I(V) characteristic occur. Addi- 
tionally, the differential conductivity at the current onset just beyond the Coulomb gap is 
calculated. 

Figure [I] shows a system of iV tunnel junctions coupled in series. In the semiclassi- 
cal model, the state of each tunnel junction is characterized by the voltage drop across 
[Fig. [TJ (a)]. The individual junction voltages can be calculated using Kirchhoff's law to- 
gether with Gauss's law. Due to tunneling there might be additional electrons on the center 
electrodes. The voltage V%. across the k-th junction with rij extra electrons on the j-th 
electrode and an externally applied voltage V can be written as 



V h {...n j ...,V) = ±- 



/ N ( m-1 \ \ 

' C X V -de E l/C m E n k _ x x 

m=2 \ j=l 
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N 
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(1) 



where C m denotes the capacitance of the m-th tunnel junction. Vk causes a mutual shift in 

the Fermi energies of the electrodes. 

The tunneling rates across the junctions, which can be calculated using Fermi's golden 
rule, are given by the following expression (for tunneling through the k-th junction from 
right to left): 
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/ + OO •/jr 
dE — \T{E)\ 2 D k _ x {E - E k ^)f(E - E k ^)D k (E - E k ) [1 - f(E - E k )\ , (2) 



where f(E) is the Fermi-Dirac distribution. If we consider, as usual, the density of states 
D near the Fermi level and the tunnel matrix element T to be energy-independent, that is 
D k (E — E k ) = D k and \T(E)\ 2 = | To | 2 , the tunneling rates from right to left (r k ) and reverse 
(Z fc ) [cf. Fig.|(b)] are 

1 A E*~ 

r k (...nj ...,V) = — - A * /? . (3) 

e 2 R k 1 - eyLp(-AEjr/k B T) 

and 

1 A 

l k (... nj ...,V) = ^p- ^ , (4) 

e 2 R k 1 - exp(— AE^/hsT) 

where the tunneling resistances R k are given by l/R k = (27ie 2 /h) Dl_ x D k \TQ\ 2 . The energy 

that an electron gains by tunneling through the fc-th junction, AE^ or AE^, may be calcu- 
lated by integrating the difference between the neighboring Fermi levels over the tunneling 
event: 



AEi 



| = -e J dq V k (. . . , n fc _i - q,n k + q,... 



V). (5) 



This is just the change in the electrostatic energy of the k-th capacitor. Evaluating the 
integral leads to a simple expression for the energy changes: 



l E l\^eV k( ...n,..,vf % (6) 

* ' Cld + Ci E i/c m ) 1L,k 

m=2 



The I(V) characteristics of this system can be determined numerically by Monte Carlo 
simulation. One main feature of the results which arises in the simulations is a widening of 
the Coulomb gap with increasing number of tunnel junctions (Fig.0). This can be understood 
for the zero-temperature limit of the system. Also, the differential conductivity at current 
onset (see inset of Fig.|3|) can be calculated in this limit. For T = 0, a tunneling event is 
impossible whenever a tunneling electron would lose energy, that is AE^T < or AE^ < 0, 
respectively. From any state defined by the numbers rij of extra electrons on the center 
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electrodes, the system may undergo a number of transitions by tunneling of an electron 
somewhere in the series. The energy changes associated with these transitions, AEJT and 
AE^, are directly related to the applied voltage via Eq.(||). In the simplest case of only one 
extra electron in the series (or one extra hole, i.e., the absence of one electron), these energy 
differences are negative for zero voltage, i.e., tunneling is impossible. An increasing voltage 
causes the energies to increase towards and then above zero, thus increasing the number 
of possible tunneling events. The energy differences associated with tunneling through the 
terminating junctions of the chain are the last to become positive. Therefore, the Coulomb 
gap of multiple junction systems can be determined by examining the voltages at which 
either AE^* and AE^~ or AE^ and AE^ cross zero. For C\ < Cjv, the first tunnel junction 
is the one that dominates the process and the Coulomb gap is found to be 

N l 

A^a P = e£- (7) 

i=2 U i 

For C\ > Cjy, the last junction dominates the process, which means that the sum in Eq.([7|) 
runs from i — 1 to N — 1. A similar approach to calculate the size of the Coulomb gap was 
chosen before, |]] although the system of junctions is not as general as that considered here 
(the capacitances are taken to be all equal). 

Higher threshold voltages can in principle be determined in a similar manner, thus ac- 
counting for the structure of multiple junction I(V) curves. Using a fairly simple computer 
program, one can calculate the states that are accessible from the "ground state" , which is 
given by nj = for all j's, and one can determine the threshold voltages at which new states 
become accessible. A state is accessible if a sequence of tunneling events with nonvanish- 
ing probability leads to it. These threshold voltages provide an intuitive understanding of 
the steps in the I(V) characteristic as well as an exact prediction of the voltages at which 
they occur. Until now, the step structures in such simulations have usually been presented 



without explanation or simply been characterized as unusual. (10]|TT| In the present analysis, 
however, the number of predicted threshold voltages exceeds the number of steps that occur 
in the numerical simulation. The occurence of steps depends crucially on the choice of the 
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resistance values, as will be discussed subsequently. 

Just beyond the Coulomb gap, i.e., at voltages slightly above the first threshold voltage 

Vtb = ^AV^ap, the most probable process is one electron (hole) tunneling through the whole 

chain of junctions before the next electron (hole) enters. The average time that the electron 

spends on one particle of the chain is given by the inverse tunneling rate. The average time 

r to tunnel through the whole chain is therefore 

1 1 1 

+ 1 7n ^7V + --- + 77T^ 7T^- ( 8 ) 



i N {o,...,o,v) i N ^(o,...,o,i,v) '■■ h(i,o,...,o,vy 

This time r is valid for voltages V slightly above the threshold voltage: V = V t h_ + 5V, where 
SV is small. At zero temperature, the tunneling rates 4 are given by the expresssion 

for AEI* < 0, 

/*(... n j ...,V) = { (9) 
[AE^/(e 2 R k ) forA£r>0. 

Therefore, just above the threshold voltage V^, the voltage at which the energy change at 
the iVth junction, 

AEJ = ^ SV, 

m=l 

crosses zero, the first term in Eq.(Q) dominates. The I(V) dependence near V t h is thus given 
by 

I(V) = - = el N (0,...,0,V) 

T 

for V < V th , 

(10) 

(C N R N E l/CJ-'dV iorV>V th . 

m=l 

In the inset of Fig.^, Monte Carlo results obtained close to the first threshold voltage are 
shown for two sets of parameters and are compared with the asymptotes given by Eq. (|T0|) . 
This figure demonstrates how drastically the differential conductivities at current onset, 
given the same average slope of 1/ Y^iLi Ri, depend on the value of the resistance of junction 
1, which in this case determines the threshold voltage. It also shows the good agreement of 
the numerical data with the analytical result. 
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In a more rigorous approach, probabilities P(. . .rij . . .) must be assgined to the states 
(. . .rij . . .) of the system. These probabilities must be solutions of the stationary master 
equation of the system. The correct expansion for the stationary current is then given by 

I = e Yl P(...n j ...)[l N (...n j ...,V)-r N (...n j ...,V)\, (11) 

all [...rij...) 

summing over all possible electronic states (. . . rij . . .) of the system with an arbitrary number 
of extra electrons on the center electrodes. 

For voltages V = V t h + 8V near the first threshold voltage, the probabilty of the "ground 
state" (0, . . . , 0) is almost 1, i.e., P(0, . . . , 0) ~ 1, and the probabilities of all other states are 
negligible. This results from the fact that, while the transition rates from the state (0, . . . , 0) 
to neighboring states are zero or nearly zero, the transition rates back to the ground state 
are already finite. Furthermore, vn(0, . . . ,0,V) vanishes for positive voltages V, leaving 
only Zjv(0, . . . , 0, V). Therefore Eq.([Tl|) reduces to 

I = el N (Q,...,Q,V) (12) 

near the first threshold voltage, which is the same result as obtained above in Eq. flTU|) . 

It is, of course, quite difficult to calculate analytically the probabilties of the states from 
the full master equation of the system, although it is possible (for capacitances all equal) 
to write an analytical expression for the average current without an explicit solution of 
the master equation in the case when the voltage corresponds to the first two steps in the 
/(V) characteristic. [0 Qualitatively, one can expect — as at the first threshold voltage 
- the steepest steps in I(V) to occur when tunneling through the junction with the lowest 
resistance value creates access to a new state. Higher steps may only occur when the 
associated threshold voltage is determined by a junction k for which 

N -l 1 

(Cfco R k E V^m) > 77 • 

m=1 E Ri 

i=l 

This tendency can also be seen in the results of Fig.^|. The solid arrows indicate voltages 
at which new states become accessible by tunneling through junction 1, which manifests itself 
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as steps in the (thicker) I(V) curve for parameter set 1 (Ri small). At voltages indicated 
by the dashed arrows, new states become accessible by tunneling through the last junction 
6, corresponding to steps in the (thinner) curve for parameter set 2 (R 6 small). 

In conclusion, Monte Carlo simulations of I(V) characteristics in a one-dimensional ar- 
rangement of tunnel junctions as well as an analytical approach have been presented. The 
analytical results may be used to understand the features found in the simulations. The re- 
sults are in general useful for the qualitative and quantitative understanding of experimental 
/(V) curves. The voltages at which steps occur can be predicted using the presented formal- 
ism — in particular, the size of the Coulomb gap can be determined. Efforts are currently 
in progress to extend the theoretical approach to finite temperatures and to investigate in 
more detail the influence of certain system parameters, such as the resistance values in the 
junction array, on the resulting I(V) curves. 
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FIGURES 

FIG. 1. (a) Circuit diagram of a system of N tunnel junctions in series, represented by their 
capacitances Cj. The voltage drop across the i-th tunnel junction is Vj. (b) Schematic diagram of 
the tunnel junctions coupled in series, showing the various tunneling rates k, r{. 

FIG. 2. Monte Carlo simulations of I(V) curves for (a) 3, (b) 5 and (c) 10 tunnel junctions. 
The parameters for the simulations are C\ = 1.2 aF, C 2 = ■ ■ ■ = Cjv-i = 1-4 aF, Cn = 2.8 aF and 
Ri = 0.01 Mfi, R 2 = --- = R N _ r = 280MO, R N = 2.8 MO and T = 4.2K. 

FIG. 3. Monte Carlo simulation and analytical predictions for the /(V) step structure of a 
6-junction system. The parameters are T = 0.01K, C\ = 1.2 aF, C2...5 = 2.8 aF, C\ = 4.4 aF (in 
both cases), R x = 24MO, R 2 ... 5 = 560MO, R 6 = 8000Mft (thick line, set 1) and i?i = 8000Mfi, 
R2...5 = 560MJ7, i?6 = 24Mf2 (thin line, set 2). Solid arrows indicate access through junction 
1, dashed arrows through junction 6 (see text). In the inset, current onsets for the two sets of 
parameters are magnified. The diamonds (set 1) and squares (set 2) are Monte Carlo results, 
whereas the solid lines display the asymptotes for V — > l 7 ^ from the analytical expression. The 
dot-dashed line shows the average slope for reference. 



Fig.l: U.E.Volmar et al., Electronic transport in a series 
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Fig.2: U.E.Volmar et al., Electronic transport in a series 




Fig.3: U.E.Volmar et al., Electronic transport in a series 
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